**********************************************************************************************************
* This replication file estimates an ordered probit model of notary entry based on the European dataset. *
**********************************************************************************************************
*1. Sample, data management and summary statistics
*2. Ordered probit entry models by country, replicating Tables 2, A3 and A4. 
*3. Calculation of notaries per capita, replicating Figure 1.

*****************************************
* 0 Set path to files and load the data *
*****************************************
clear
clear matrix
clear mata
set more off
cd $main_directory
use "Data_Europe\Generated_data\country_data.dta", clear
drop if Population==0

bysort Country: su Population Area_sqkm

*****************************************************
* 1. Sample, data management and summary statistics *
*****************************************************


*Select sample
gen Density_th=Population/(1000*Area_sqkm)
drop if Density_th>1.2

*Data management
replace NotOffices=4 if NotOffices>4

gen lnPop=log(Population)
gen lnArea=log(Area_sqkm)

gen IncMun_th=INCOME/1000

bysort Country: egen temp=mean(IncMun_th)
replace IncMun_th=temp if IncMun_th==0

bysort Country: egen m_Population=mean(Population)
bysort Country: egen m_IncMun_th=mean(IncMun_th)
bysort Country: egen m_Perc_Male=mean(Perc_Male)
bysort Country: egen m_Perc_Young=mean(Perc_Young)
bysort Country: egen m_Perc_Elderly=mean(Perc_Elderly)
bysort Country: egen m_lnArea=mean(lnArea)

*Summary statistics: Table A3, Panel A
sum lnPop IncMun_th Perc_Male Perc_Young Perc_Elderly lnArea, d 
*Summary statistics: Table A3, Panel B
tab Country NotOffices


*********************************************
* 2. Ordered probit entry models by country *
*********************************************

global xvars "lnPop lnArea IncMun_th Perc_Male Perc_Young Perc_Elderly "
label variable lnPop "Population (log)"
label variable IncMun_th "Mean income in 1000 \euro{}"
label variable Perc_Male "\% male inhabitants"
label variable Perc_Young "\% young inhabitants"
label variable Perc_Elderly "\% elderly inhabitants"
label variable lnArea "Area (log)"

* start with AT to initialize variables, then do rest in foreach loop
qui gen selected_country="AT"
qui oprobit NotOffices $xvars if Country==selected_country
est sto oprob_AT
qui gen et01_notaries=round(exp(-(_b[IncMun_th]*m_IncMun_th+_b[Perc_Male]*m_Perc_Male+_b[Perc_Young]*m_Perc_Young+_b[Perc_Elderly]*m_Perc_Elderly+_b[lnArea]*m_lnArea-_b[/cut1])/_b[lnPop])) if Country==selected_country
qui gen et12_notaries=round(exp(-(_b[IncMun_th]*m_IncMun_th+_b[Perc_Male]*m_Perc_Male+_b[Perc_Young]*m_Perc_Young+_b[Perc_Elderly]*m_Perc_Elderly+_b[lnArea]*m_lnArea-_b[/cut2])/_b[lnPop])/2) if Country==selected_country
qui gen et23_notaries=round(exp(-(_b[IncMun_th]*m_IncMun_th+_b[Perc_Male]*m_Perc_Male+_b[Perc_Young]*m_Perc_Young+_b[Perc_Elderly]*m_Perc_Elderly+_b[lnArea]*m_lnArea-_b[/cut3])/_b[lnPop])/3) if Country==selected_country
qui gen et34_notaries=round(exp(-(_b[IncMun_th]*m_IncMun_th+_b[Perc_Male]*m_Perc_Male+_b[Perc_Young]*m_Perc_Young+_b[Perc_Elderly]*m_Perc_Elderly+_b[lnArea]*m_lnArea-_b[/cut4])/_b[lnPop])/4) if Country==selected_country
nlcom (etr12_not:exp((_b[/cut2]-_b[/cut1])/_b[lnPop])*1/2) (etr23_not:exp((_b[/cut3]-_b[/cut2])/_b[lnPop])*2/3) (etr34_not:exp((_b[/cut4]-_b[/cut3])/_b[lnPop])*3/4),post

qui gen etr12_notaries=round(_b[etr12_not],0.001) if Country==selected_country
qui gen etr12_notaries_se=round(_se[etr12_not],0.001) if Country==selected_country
qui gen etr23_notaries=round(_b[etr23_not],0.001) if Country==selected_country
qui gen etr23_notaries_se=round(_se[etr23_not],0.001) if Country==selected_country
qui gen etr34_notaries=round(_b[etr34_not],0.001) if Country==selected_country
qui gen etr34_notaries_se=round(_se[etr34_not],0.001) if Country==selected_country

qui test _b[etr12_not]=1
qui gen etr12_notaries_p_val=round(r(p),0.001) if Country==selected_country
qui test _b[etr23_not]=1
qui gen etr23_notaries_p_val=round(r(p),0.001) if Country==selected_country
qui test _b[etr34_not]=1
qui gen etr34_notaries_p_val=round(r(p),0.001) if Country==selected_country
est sto pr_AT

local countrynames "BE DE FR IT NL PT ES"
foreach x of local countrynames {
qui replace selected_country="`x'"
qui oprobit NotOffices $xvars if Country==selected_country
est sto oprob_`x'
qui replace et01_notaries=round(exp(-(_b[IncMun_th]*m_IncMun_th+_b[Perc_Male]*m_Perc_Male+_b[Perc_Young]*m_Perc_Young+_b[Perc_Elderly]*m_Perc_Elderly+_b[lnArea]*m_lnArea-_b[/cut1])/_b[lnPop])) if Country==selected_country
qui replace et12_notaries=round(exp(-(_b[IncMun_th]*m_IncMun_th+_b[Perc_Male]*m_Perc_Male+_b[Perc_Young]*m_Perc_Young+_b[Perc_Elderly]*m_Perc_Elderly+_b[lnArea]*m_lnArea-_b[/cut2])/_b[lnPop])/2) if Country==selected_country
qui replace et23_notaries=round(exp(-(_b[IncMun_th]*m_IncMun_th+_b[Perc_Male]*m_Perc_Male+_b[Perc_Young]*m_Perc_Young+_b[Perc_Elderly]*m_Perc_Elderly+_b[lnArea]*m_lnArea-_b[/cut3])/_b[lnPop])/3) if Country==selected_country
qui replace et34_notaries=round(exp(-(_b[IncMun_th]*m_IncMun_th+_b[Perc_Male]*m_Perc_Male+_b[Perc_Young]*m_Perc_Young+_b[Perc_Elderly]*m_Perc_Elderly+_b[lnArea]*m_lnArea-_b[/cut4])/_b[lnPop])/4) if Country==selected_country
qui nlcom (etr12_not:exp((_b[/cut2]-_b[/cut1])/_b[lnPop])*1/2) (etr23_not:exp((_b[/cut3]-_b[/cut2])/_b[lnPop])*2/3) (etr34_not:exp((_b[/cut4]-_b[/cut3])/_b[lnPop])*3/4),post

qui replace etr12_notaries=round(_b[etr12_not],0.001) if Country==selected_country
qui replace etr12_notaries_se=round(_se[etr12_not],0.001) if Country==selected_country
qui replace etr23_notaries=round(_b[etr23_not],0.001) if Country==selected_country
qui replace etr23_notaries_se=round(_se[etr23_not],0.001) if Country==selected_country
qui replace etr34_notaries=round(_b[etr34_not],0.001) if Country==selected_country
qui replace etr34_notaries_se=round(_se[etr34_not],0.001) if Country==selected_country

qui test _b[etr12_not]=1
qui replace etr12_notaries_p_val=round(r(p),0.001) if Country==selected_country
qui test _b[etr23_not]=1
qui replace etr23_notaries_p_val=round(r(p),0.001) if Country==selected_country
qui test _b[etr34_not]=1
qui replace etr34_notaries_p_val=round(r(p),0.001) if Country==selected_country
est sto pr_`x'
}


*Ordered probit results: Table A4
estimates replay oprob_AT oprob_BE oprob_FR oprob_DE oprob_IT oprob_NL oprob_PT oprob_ES
esttab oprob_AT oprob_BE oprob_FR oprob_DE oprob_IT oprob_NL oprob_PT oprob_ES using "Estimates/Table_A4.tex", replace ///
		b(%8.3f) se(%8.3f) scalars(N) sfmt(0) star(* 0.10 ** 0.05  *** 0.01) label compress nodepvar nonotes ///
		nogap collabels(none) eqlabels(none)  ///
		title(Estimation results from ordered probit entry model \label{tab:entry_europe})     ///
		nonumbers mtitles( "Austria" "Belgium" "France" "Germany" "Italy" "Netherlands" "Portugal" "Spain") 

*Entry thresholds and ratios: Table 2
estimates replay pr_AT pr_BE pr_DE pr_FR pr_IT pr_NL pr_PT pr_ES
preserve

collapse  et01_notaries et12_notaries et23_notaries et34_notaries etr12_notaries etr12_notaries_p_val etr23_notaries etr23_notaries_p_val etr34_notaries etr34_notaries_p_val, by(Country)
dataout, save("Estimates/Table_2.tex") tex replace
restore

*****************************************
* 3. Calculation of notaries per capita *
*****************************************
use "Data_Europe\Generated_data\country_data.dta", clear


preserve
collapse (sum) Population NotOffices Notaries, by(Country)
gen order=.
replace order=1 if Country=="BE"
replace order=2 if Country=="DE"
replace order=3 if Country=="IT"
replace order=4 if Country=="FR"
replace order=5 if Country=="NL"
replace order=6 if Country=="AT"
replace order=7 if Country=="ES"
replace order=8 if Country=="PT"
gen NotOff_Inhab=NotOffices/(Population/100000)
gen Not_Inhab=Notaries/(Population/100000)
//Figure 1
graph bar NotOff_Inhab Not_Inhab, over(Country, sort(order)) bar(1, color(gs5)) bar(2, color(gs14)) ///
    graphregion(color(white)) bgcolor(white) ///
    ytitle("Number of sellers per 1000 inhabitants") ///
    legend(label(1 "Notary offices per capita") label(2 "Notaries per capita"))
	graph export ".\Estimates\Graphs\Figure1.png", replace
	graph export "Estimates\Graphs\Figure_1.eps", replace

restore
